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Abstract — We  present  a  novel  two-stage  methodology  for 
locating  a  Chemical,  Biological,  Radiological,  or  Nuclear 
(CBRN)  source  in  an  urban  area  using  a  network  of 
sensors.  In  contrast  to  earlier  work,  our  approach  does 
not  solve  an  inverse  dispersion  problem  but  relies  on  data 
obtained  from  a  simulation  of  the  CBRN  dispersion  to 
obtain  probabilistic  descriptors  of  sensor  measurements 
under  a  variety  of  CBRN  release  scenarios.  At  its  first  stage, 
subsequent  sensor  observations  under  nominal,  CBRN 
event-free  conditions  are  assumed  to  be  independent  and 
identically  distributed  and  we  rely  on  the  method  of  types 
to  detect  a  CBRN  event.  Conditional  on  such  an  event, 
subsequent  sensor  observations  are  assumed  to  follow  a 
Markov  process.  Using  composite  hypothesis  testing  we 
map  sensor  measurements  to  a  source  location  chosen  out 
of  a  discrete  set  of  possible  locations.  We  leverage  large 
deviation  techniques  to  obtain  a  bound  on  the  localization 
probability  of  error  and  propose  several  methodologies 
for  fusing  sensor  data  to  arrive  at  a  localization  decision, 
including  a  distributed  one.  We  also  address  the  problem 
of  optimally  placing  sensors  to  minimize  the  localization 
probability  of  error.  Our  techniques  are  validated  numer¬ 
ically  using  two  different  CBRN  release  simulators. 

Index  Terms — Source  detection,  source  localization,  com¬ 
posite  hypothesis  testing,  large  deviations,  optimization, 
sensor  placement. 

I.  Introduction 

CONCERN  regarding  Chemical,  Biological,  Radio¬ 
logical,  or  Nuclear  (CBRN)  terrorism  is  steadily 
increasing.  This  is  due  to  gains  in  the  technological 
capabilities  of  producing  existing  and  potentially  new, 
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more  lethal  CBRN  agents  and  delivery  mechanisms. 
Urban  areas  are  of  particular  concern  since  these  areas 
tend  to  have  large  population  densities  and  are  centers 
for  large-scale  commerce  and  politics  [1], 

In  this  context,  it  becomes  critical  to  detect  and  locate 
a  CBRN  source.  Early  detection  and  accurate  localiza¬ 
tion  enable  effective  emergency  response.  Localization 
is  important  even  for  weak  CBRN  sources  as  it  can  help 
prevent  a  harmful  release  and  facilitate  investigating  the 
origins  of  the  source  (forensics). 

Invariably,  all  existing  approaches  rely  on  a  network  of 
sensors  deployed  in  the  urban  area  and  able  to  measure 
CBRN  particulate  concentrations.  They  differ  however, 
in  how  they  process  the  sensor  data.  A  related  problem  is 
that  of  placing  a  given  set  of  sensors  so  as  to  optimize  a 
performance  metric  related  to  detection  and  localization. 

Earlier  work  (reviewed  later  in  more  detail)  has  mostly 
concentrated  on  solving  a  difficult  inverse  dispersion 
problem  in  order  to  locate  the  source.  Such  an  approach 
requires  an  explicit  model  of  CBRN  dispersion  and, 
typically,  simple  analytical  models  are  being  assumed 
to  render  the  inverse  problem  tractable.  Yet,  simple 
dispersion  models  are  not  accurate  in  urban  areas  for 
a  variety  of  reasons.  First,  a  typical  city  tends  to  be 
characterized  by  irregular  geometry,  many  different  types 
of  structures  and  buildings  (each  with  its  own  surface 
texture),  and  highly  dynamic  population  fluxes,  thus, 
giving  rise  to  an  inhomogeneous  dispersion  terrain.  Sec¬ 
ond,  weather  patterns,  which  drive  CBRN  dispersion, 
are  highly  variable  and  interact  with  the  various  city 
structures  to  create  complex  phenomena  (micro-climate 
effects,  urban  canyons,  etc.).  All  this  complexity  can  not 
be  captured  analytically  and,  coupled  with  the  inherent 
uncertainty,  suggests  a  stochastic  approach. 

In  this  paper  we  develop  such  an  approach  that  uses 
simulation  of  CBRN  dispersion  under  a  variety  of  source 
locations  and  dispersion  scenarios  to  develop  probabilis¬ 
tic  descriptors  of  sensor  measurements.  Localization  is 
performed  solely  through  consultation  of  the  sensor  mea¬ 
surements,  not  through  an  explicit  physical  dispersion 
model.  Hence,  we  describe  our  approach  as  model-free. 
Our  approach  consists  of  two  stages:  (a)  detection,  and 
(6)  localization.  To  detect  a  source  we  assume  that  subse¬ 
quent  measurements  at  any  given  sensor  are  independent 
and  identically  distributed  (iid).  This  is  justified  since,  in 
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the  absence  of  a  source,  sensor  measurements  are  driven 
by  noise.  Using  the  method  of  types  ([2],  [3])  we  com¬ 
pare  sequences  of  measurements  with  the  probabilistic 
descriptors  derived  from  the  simulation  to  determine  the 
presence  of  a  source.  We  develop  a  rigorous  detection 
test  and  show  asymptotic  Neyman-Pearson  optimality. 

Positive  source  detection  triggers  the  2nd  -  localiza¬ 
tion  -  stage  of  our  approach.  During  that  stage,  and 
conditional  on  the  presence  of  particulates  from  a  source, 
subsequent  measurements  at  a  sensor  are  dependent. 
We  make  a  Markovian  assumption  and  develop  an  ap¬ 
proach  to  compare  measurements  from  the  sensors  to  the 
probabilistic  descriptors  obtained  from  the  simulation. 
Localization  is  formulated  as  a  composite  hypothesis 
testing  problem.  We  devise  several  versions  of  a  local¬ 
ization  algorithm  -  a  distributed,  a  centralized,  and  a 
hybrid  -  using  different  forms  of  data  fusion  from  the 
multiple  sensors.  We  derive  a  bound  on  the  probability 
of  error  by  extending  our  earlier  work  from  [4],  which 
also  considered  a  localization  problem  but  in  a  very 
different  application  domain.  Further,  we  use  this  bound 
to  formulate  the  problem  of  placing  sensors  so  as  to  min¬ 
imize  the  probability  of  error.  As  our  numerical  results 
demonstrate  such  an  optimization  can  lead  to  a  dramatic 
improvement  in  the  error  probability  compared  to  ad-hoc 
sensor  placement.  The  end  result  of  our  approach  is  a 
placement  of  the  sensors  and  an  asymptotic  guarantee 
on  the  corresponding  probability  of  error. 

One  will  notice  some  similarities  between  the  theoret¬ 
ical  material  presented  in  [4]  and  in  this  paper.  We  recast 
theorems  and  definitions  here  in  a  new  application  both 
for  the  sake  of  completeness  as  well  as  to  provide  back¬ 
ground  for  our  (Markovian)  extensions  to  the  previous 
(iid)  work. 

We  summarize  our  contributions  below: 

1)  We  propose  a  new  paradigm  for  CBRN  source  de¬ 
tection  and  localization  that  bypasses  many  of  the 
problems  existing  “inverse  dispersion”  approaches 
face  in  an  inhomogeneous  and  uncertain  urban 
setting. 

2)  We  extend  our  earlier  asymptotic  performance 
analysis  of  composite  hypothesis  testing  in  [4] 
from  the  iid  to  a  Markovian  setting. 

3)  We  demonstrate  how  to  place  sensors  so  as  to 
minimize  the  localization  probability  of  error. 

4)  We  validate  our  approach  using  two  different, 
and  independently  developed,  dispersion  simula¬ 
tion  engines:  the  Quick  Urban  &  Industrial  Com¬ 
plex  (QUIC)  Dispersion  Modeling  System  [5]  and 
a  new  simulator  based  on  the  Lattice  Boltzmann 
method  ( LBM )  that  we  develop  for  our  purposes. 
The  latter  simulator  is  of  independent  interest 
as  it  is  known  to  be  a  naturally  parallelizable 
algorithm  [6],  This  quality  places  large  scenarios 


computationally  within  reach.  These  simulators 
simply  simulate  the  physics  of  dispersion  from 
some  point  release  under  given  weather  conditions 
and  make  none  of  the  distributional  assumptions 
our  approach  posits.  Although  it  is  certainly  true 
that  QUIC  and  LBM  make  certain  assumptions 
and  approximations  to  simplify  the  physics  of  dis¬ 
persion,  our  approach  decouples  these  assumptions 
from  the  process  of  release  detection  and  source 
localization  by  relying  solely  on  measurements  ob¬ 
tained  from  sensors  during  particulate  propagation. 
As  such,  our  approach  naturally  extends  to  any 
simulator  or  even  to  measurements  obtained  from 
physical  experimentation.  Different  scenarios  will 
call  for  different  dispersion  simulators.  If  accuracy 
is  less  of  an  issue,  one  could  select  a  faster,  less 
accurate  simulator.  If  the  area  under  observation 
contains  heavy  foliage,  a  simulator  with  satis¬ 
factory  deposition  accuracy  should  be  selected. 
For  situations  in  which  particulate  deposition  is 
important  due  to  physicochemical  properties  of 
the  dispersed  agent  or  sensor  performance  limi¬ 
tations,  say,  the  computationally  efficient  method 
of  estimating  boundary  properties  presented  in  [7] 
could  be  of  use.  A  set  of  illustrative  numerical 
results  demonstrate  excellent  performance  of  our 
methods.  The  results  show  promise  for  a  practical 
CBRN  urban  detection  and  localization  system 
based  on  this  work. 

The  remainder  of  the  paper  is  organized  as  follows. 
We  discuss  related  work  in  Section  II.  We  formulate 
our  problem  in  Section  III  where  we  also  introduce 
some  of  our  notation.  We  deal  with  the  detection  stage 
in  Section  IV  and  the  localization  stage  in  Section  V. 
The  placement  problem  is  discussed  in  Section  VI. 
The  two  dispersion  simulators  and  a  host  of  numerical 
results  are  presented  in  Section  VII.  Final  remarks  are 
in  Section  VIII.  The  Appendix  contains  technical  proofs 
of  our  results. 

Notational  Conventions:  Throughout  the  paper  all 
vectors  are  assumed  to  be  column  vectors.  We  use  lower 
case  boldface  letters  to  denote  vectors  and  for  economy 
of  space  we  write  x  =  (aq ,...,xr)  for  the  column 
vector  x.  x'  denotes  the  transpose  of  x  and  0  the 
vector  of  all  zeroes.  We  use  upper  case  boldface  letters 
to  denote  matrices.  Finally,  we  use  \g/\  to  denote  the 
cardinality  of  a  set  g/. 

II.  Related  work 

A  sequential  measurement  based  hypothesis  testing 
paradigm  for  CBRN  release  detection  is  presented  in  [8], 
There,  the  authors  make  the  assumption  that  concentra¬ 
tion  samples  adhere  to  a  zero-mean  Gaussian  error  and 
formulate  a  maximum  likelihood  test  used  for  release 
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detection.  Our  detection  approach  is  similar  in  that  it  is 
a  hypothesis  test.  However,  we  make  no  assumptions 
on  the  error  a  CBRN  sensor  makes,  opting  instead 
to  learn  a  probabilistic  law  of  nominal,  release-free 
sensor  observations.  Whenever  a  sensor’s  observations 
sufficiently  depart,  according  to  a  Hoeffding  test,  from 
the  known  nominal  behavior,  we  declare  that  a  CBRN 
release  has  occurred. 

Existing  CBRN  source  localization  approaches  ([9], 
[10],  [11])  observe  CBRN  agent  concentrations  and  solve 
the  inverse  problem  of  tracing  dispersion  backward  in 
time  and  space  to  the  source  of  the  release.  As  discussed, 
limitations  of  this  methodology  stem  from  the  irregular 
and  dynamic  phenomena  typically  found  in  urban  areas. 
In  [9]  the  presence  of  challenging  geographies  and  wind 
turbulence  is  accommodated  by  incorporating  Monte 
Carlo  simulation  of  fluid  dispersion.  However,  all  of 
these  works  suffer  from  the  difficulty  of  determining, 
without  a  detection  process,  the  point  in  time  in  which 
the  event  started.  Barring  this  information,  these  inverse 
problem  approaches  are  vulnerable  to  erroneous  local¬ 
izations. 

Simulation-aided  localization  of  contaminant  releases 
on  a  continental  scale  is  presented  in  [12].  In  [13],  the 
inverse  problem  is  replaced  by  computing  the  posterior 
probabilities  of  a  finite  set  of  release  locations  and 
selecting  the  release  location  with  maximal  likelihood. 
The  computed  likelihoods  are  based  on  assumptions 
on  both  the  distribution  of  concentration  observations 
at  sensor  locations  and  a  particular  dispersion  model. 
Our  approach  is  similar  to  these,  in  that  localization  is 
based  entirely  on  sensor  observations  and  probabilistic 
descriptors  of  CBRN  releases  obtained  from  simulation, 
with  the  important  difference,  however,  that  the  only 
assumptions  we  make  on  concentration  observations 
under  release  conditions  is  that  they  adhere  to  a  first- 
order  Markov  chain. 

While  our  approach  does  depend  on  numerical  or 
physical  simulation  of  agent  dispersion  to  build  a  sta¬ 
tistical  understanding  of  CBRN  release  evolution,  once 
this  understanding  is  established  the  actual  dispersion 
model  is  no  longer  used.  Hence,  we  refer  to  our  lo¬ 
calization  approach  as  model-free.  Additionally,  while 
we  employ  two  specific  dispersion  simulators  in  the 
numerical  evaluation  in  Section  VII,  any  dispersion 
simulation,  numerical  or  physical,  can  be  used.  We 
suggest  using  the  simulation  scheme  that  provides  the 
most  accurate  depiction  of  particulate  releases  for  the 
particular  application.  This  is  even  more  important  in 
applications  where  turbulence  plays  a  greater  role  in 
particulate  displacement. 

Our  localization  approach  is  inspired  by  [4].  There, 
a  wireless  sensor  deployed  in  an  indoor  environment  is 
located  through  sequential  hypothesis  testing  by  using 


signal  strength  measurements  corresponding  to  packets 
transmitted  by  the  sensor  and  received  by  a  set  of  station¬ 
ary  clusterheads.  [4]  also  solves  the  problem  of  optimally 
placing  clusterheads  to  reduce  erroneous  localizations. 
The  key  difference  between  that  work  and  the  work 
presented  here,  aside  from  vastly  different  application 
areas,  is  that  subsequent  observations  made  by  the  clus¬ 
terheads  in  the  wireless  sensor  localization  problem  can 
be  considered  iid.  Under  a  CBRN  release  scenario,  if  a 
CBRN  sensor  observes  large  particle  concentrations  it  is 
likely  to  observe  large  particle  concentrations  in  its  next 
sample.  To  that  end,  in  this  work  we  extend  the  theory 
presented  in  [4]  to  allow  for  dependencies  in  subsequent 
sensor  observations. 

Our  optimal  sensor  placement  approach  also  stems 
from  the  work  presented  in  [4],  which  poses  clusterhead 
placement  as  an  optimization  problem  whose  objective 
is  to  minimize  the  probability  of  localization  error. 
This  differs  from  the  more  traditional  sensor  placement 
methodologies  [14],  [15],  whose  objectives  are  the  facil¬ 
itation  and  maintenance  of  general  observability  under 
the  possibility  of  sensor  faults.  To  our  knowledge  there 
are  no  other  existing  CBRN  sensor  placement  techniques 
that  are  motivated  by  analytically  derived  localization 
error  probabilities. 

III.  Problem  formulation 

Consider  a  CBRN  sensor  network  deployed  in  an 
urban  setting.  Due  to  the  nature  of  plume  dispersion 
and  existing  CBRN  attack  response  techniques,  extreme 
precision  on  the  locations  of  releases  and  sensors  is 
not  needed.  Rather,  general  locations  (e.g.,  corner  of 
X  St.  and  Y  Ave.)  suffice  for  disaster  avoidance  and 
response  measures.  It  naturally  follows  that  release  and 
potential  sensor  locations  can  be  discretized  to  conform 
to  potentially  irregular  grids.  In  the  following,  we  assume 
N  possible  CBRN  release  scenarios  represented  by  the 
set  2z?  =  {Li, . . . ,  Ljv}  and  M  possible  CBRN  sensor 
locations  represented  by  the  set  3$  =  { B\ , . . . .  Bm  }. 
For  applications  in  which  greater  precision  is  required, 
the  set  of  discrete  locations  can  be  extended  to  include 
more  locations.  However,  this  accrues  the  cost  of  requir¬ 
ing  more  simulations.  The  elements  of  Jz?  correspond  to 
a  release  location  and  a  set  of  release  associated  charac¬ 
teristics  potentially  including  wind  direction,  variability 
in  wind  speed,  and  release  mass,  among  others.  Wind,  in 
particular,  is  a  critical  component  of  a  release  scenario 
as  it  is  the  primary  mode  of  particulate  transport.  For 
a  release  scenario  L;/ ,  the  corresponding  set  of  release 
characteristics  will  be  denoted  f h  =  {6], . . .  ,6^^}. 

Fet  y  denote  the  vector  of  sensor  observations  made 
by  a  sensor  at  location  Bk .  Each  sensor  outputs  a  CBRN 
agent  concentration  estimation  and  potentially  local  wind 
observations.  Since  different  concentrations  result  in 
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different  levels  of  casualties,  counter  measures,  and 
responses,  observed  concentrations  can  naturally  be  dis¬ 
cretized  to  concentration  levels  (e.g.,  [0,  x],  [x,ulct5o\, 
etc.,  where  t/lctso  represents  a  concentration  that  would 
lead  to  an  accumulated  dosage  that  produces  a  50% 
chance  of  survival).  Each  yk  is  mapped  to  a  symbol 
in  a  finite  alphabet  denoted  by  E  =  {<7i, . . . ,  criEi }. 
A  sequence  of  n  sensor  observations  at  location  If.  is 
denoted  yfc,n  =  (y^, . . .  ,yk).  It  is  assumed  throughout 
the  sequel  that  an  initial  state  of  y§  is  known  a  priori. 

While  not  under  an  attack  scenario,  a  sensor’s  agent 
concentration  approximation  fluctuations  are  attributed 
to  measurement  noise.  Accordingly,  when  there  is  no 
CBRN  attack,  the  sensor’s  readings  can  be  assumed  to 
be  iid. 

Under  a  release  scenario,  however,  the  agent  disperses 
gradually,  resulting  in  concentration  observations  that 
change  much  slower  than  a  sensor’s  sampling  rate.  That 
is,  if  a  certain  concentration  level  is  observed  by  a 
sensor  at  a  certain  time  step,  one  can  expect  to  see  a 
similar  concentration  observation  in  the  next  time  step. 
To  capture  this  dependency,  we  model  concentration 
level  at  a  sensor  when  there  is  an  attack  as  a  first-order 
Markov  chain. 

The  first  problem  we  wish  to  solve  is  to  detect  a 
CBRN  event,  given  a  set  of  observations  {y1’",..., 
yK'n}  from  K  sensors  placed  at  locations  B\, ... ,  Bk- 
Once  a  CBRN  event  is  identified,  the  next  problem  is  to 
determine  where  the  event  originated.  We  also  wish  to 
solve  the  related  problem  of  selecting  the  locations  from 
!■'/>  where  we  can  place  K  sensors  (with  K  <  M)  such 
that  we  minimize  the  probability  of  localization  error. 
The  solution  to  this  last  problem  is  accomplished  through 
analysis  of  theoretical  performance  guarantees. 

IV.  Release  detection 

In  this  section  we  describe  our  approach  to  the  prob¬ 
lem  of  release  detection.  We  first  construct  empirical 
probability  laws  for  the  vector  yk  under  nominal  con¬ 
ditions  and  then  compare  these  laws  with  empirical 
probability  measures  based  on  in  situ  sensor  observa¬ 
tions.  When  a  measure  differs  significantly  from  the 
established  probability  laws,  a  release  is  declared  and 
the  localization  engine  is  triggered. 

Allow  Mi(E)  to  denote  the  set  of  all  probabil¬ 
ity  vectors  with  support  defined  by  the  state  alphabet 
E  =  {or, . . . ,  CT|s|}.  Consider  the  empirical  probability 
measure  of  a  sequence  yk,n  =  (yf, . . . ,  yk), 

£yk  n  (<7r)  =  -  =  °>},  r  =  1, . . . ,  |E|,  (1) 

k  n 

where  1  { • }  denotes  the  indicator  function.  Let  <%y  = 

[Sy  ’  (or), . . . ,  Sy  (<T|ej)).  First  we  construct  the  em¬ 
pirical  probability  law  pr  £  Mi(E)  for  all  sensor 


locations  B £  SS  from  observed  sequences  yk'N  = 
(y?, . . .  ,y^)  with  sufficiently  large  N  by  setting  pk  = 
gyk,N . 

When  the  sensor  network  is  operational,  a  sequence 
yk.n  _  _  _  _  ;y js  observed  for  some  n  <C  N.  An 

empirical  probability  measure  vk  =  <§y  "  is  constructed 
via  (1).  To  determine  if  the  observations  in  yk'n  are 
anomalous,  we  employ  the  so-called  method  of  types. 
The  following  theorem  from  [3]  is  useful. 


Theorem  IV.l  For  every  v  £  Mi(E)  let 
h(y)=H{y\lxk\ 

where  H(i>\pk)  is  the  relative  entropy  of  the  probability 
vector  v  with  respect  to  pk: 


|S| 


H{v\nk)  =  ^iz(ar)log 


r=  1 


^(oy) 

pk{ar)' 


Then,  for  any  set  T  of  probability  vectors  in  M\  (E) 


—  inf  /i(iz)  <  liminf  —  logP„fc 
n—> oo  n 

1  r  fc  n 

limsup  —  logP„fc  <?y  ST 

n — »oo  Tl 

where  T°  denotes  the  interior  of  T. 


sy  £  r 


< 


<  -  inf  ii(iz), 

I'GT 


Theorem  IV.l  rigorously  identifies  a  distance  metric 
that  can  be  used  in  what  is  known  as  the  Hoeffding’s 
rule  [3]: 

{k  7i 

no  attack,  if  )  <  77,  ^ 

attack,  otherwise. 

Here,  77  is  a  threshold  that  can  be  determined  as  77  = 
—  log  - ,  where  e  is  a  tolerable  false  alarm  rate  [4] . 

In  order  to  assess  the  validity  and  optimality  of  tests  of 
the  form  (2),  we  need  to  define  what  notion  of  optimality 
we  hope  to  attain.  One  often  adopted  measure  of  any 
decision  test’s  optimality  is  to  consider  the  likelihood  of 
making  an  error.  The  two  types  of  error  the  Hoeffding 
test  in  (2)  can  make  are:  declaring  an  anomaly  when  in 
fact  the  activity  recorded  in  yk,n  was  nominal  (Type-I 
error),  or  not  declaring  an  anomaly  when  yk’n  departs 
from  nominal  behavior  (Type-II  error).  We  denote  the 
probabilities  of  these  two  errors  by  ak,n  and  /3k, n, 
respectively. 

A  particular  challenge  arises  when  trying  to  evaluate 
(3k, n ■  While  the  probability  law  depicting  nominal  behav¬ 
ior  pk  can  be  used  to  evaluate  ak,n  =  [h  (&y  " )  > 
77],  computation  of  [3k, n  requires  considering  all  proba¬ 
bility  vectors  from  the  set  Mi(E)  \  {/xfc}.  With  this  in 
mind,  we  introduce  the  following  concept  of  asymptotic 
optimality  known  as  the  Neyman-Pearson  criterion. 
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Definition  1 

A  test  SF  is  considered  optimal  for  a  given  T)  >  0  if 
among  all  tests  that  satisfy 

limsup  -otk ,n  <  ~V 

n — »oo  ^ 

it  maximizes  —  limsupn_>00  ^/3fe,n  uniformly  over  all 
probability  vectors  in  M\  (E)  \  {pk}. 

The  following  theorem,  thanks  to  Hoeffding  [16]  and 
appearing  in  [3],  establishes  the  optimality  of  (2). 

Theorem  IV.2  The  anomaly  detector  presented  in  (2)  is 
optimal  under  Definition  1. 

The  effect  of  the  first  stage  conclusion  is  two-fold. 
First,  it  establishes  a  reasonable  indication  that  a  CBRN 
event  has  actually  happened  and  that  release  localization 
analysis  is  within  order.  In  addition,  it  provides  an 
estimation  on  the  time  of  the  CBRN  release.  In  practice, 
the  time  of  attack  is  taken  to  be  the  time  at  which  the  first 
high  concentration  observation  within  y  ■  is  sampled. 
This  is  an  important  contribution  to  the  second  stage 
problem,  as  knowledge  of  the  time  of  attack  decreases 
the  amount  of  parametrizations  in  Qj  that  need  to  be 
considered.  For  localization,  the  window  of  observations 
beginning  at  the  time  of  the  release  and  of  length 
determined  by  an  ideal  response  time  is  used. 

V.  Source  localization 

Allow  M2(E  x  E)  to  denote  the  set  of  all  discrete 
time  Markov  transition  probability  matrices  on  the  state 
alphabet  E  =  {or, . . . , cr|s| }-  For  each  sensor  location 
and  release  scenario  pair  associate  a  series  of 

first-order  Markov  transition  probability  matrices  Hg., 
for  all  Oj  £  f lj,  defined  as 

n  =  {<K,^)}'S|  ,  (3) 

J  {  J  )  u,v= 1 

where  tt\.  (av,  au)  =  Pgjy tfc+1  =  au|yffc  =  av].  These 
matrices  depict  the  CBRN  particulate  concentration  state 
evolution  over  all  considered  parametrizations  of  release 
scenarios  and  available  sensor  locations.  Using  (3),  the 
probability  of  observing  a  sequence  of  sensor  observa¬ 
tions  yk’n  can  be  expressed  as 

n 

Pi  (yk,n)  =  (yk0)  I]  <  (yh,yth  (4) 

t=l 

where  Pg.(-)  denotes  the  stationary  probability  distri¬ 
bution  on  the  observed  concentration  states  at  sensor 
location  Bj.  under  release  parametrization  (f  of  scenario 
Ly 

Under  the  assumption  of  a  Markov  source,  the  em¬ 
pirical  measure  of  a  sequence  y”  =  (y0,  yi, . . . ,  y„)  of 


sensor  observations  (note  we  now  assume  that  there  is  a 
known  initial  concentration  state  yo)  takes  the  form 

n  1  n 

Sf  {av,au)  =  -  Y]  l{yi_i  =  av,yz  =  au}.  (5) 

i= 1 

For  any  release  location  and  release  scenario  pair 
( Bk,Lj )  we  can  observe  yfe,JV  under  release  conditions 
Oj  with  sufficiently  large  N  and  obtain  a  (Markovian) 
empirical  probability  law  as  follows: 

nk  =  Uf  N  {ou\<jA 

3  t  J  <Tu,<r„e£ 

where  ( cru\<jv )  are  the  empirical  transition  probabil¬ 

ities  observed  in  y”,  defined  as 


Sf  (au\av) 


Sf  (<Jv,cru) 


vn 

The  empirical  distribution  is  said  to  be  shift  invari¬ 
ant  if  the  empirical  marginal  probabilities  of  each  state 
are  equal.  That  is,  is  shift  invariant  if 


|E|  „  |s| 

(crv)  =  '^2^2  (<7v,(7u)  =  ^2^2  (cr„,crw). 


The  next  theorem,  a  version  of  Sanov’s  Theorem  for 
Markov  chains  proven  in  [3],  provides  a  convex  relative 
entropy  measure  for  Markov  chains. 


Theorem  V.l  For  every  shift  invariant  Q  £  M2( E  x  E) 
let 

|S| 

J2(Q|II)  =  ^  q(av)H2(qf\av)\ir(av,  ■)), 

V  =  1 

where  H2  is  the  relative  entropy 

|E|  (II 

^  TT(<7v,Ou) 

Then  for  any  set  F  C  M2(E  X  E)  of  shift  invariant 
Markov  sources 

-  inf  I2(Q|II*  )  <  liminf  1  logPj  Uf  "  £  T  < 
Qer°  3  n—>oo  n  3  L  J 

limsup  —  log Pg  Uf’"  erj  <  -  inf  I2(Q|nl), 
n — ^oo  n  3  L  J  Qer  3 

where  T°  denotes  the  interior  ofT  and  Pg  is  a  proba- 

i  #  3 

bility  evaluated  assuming  that  yK'n  is  drawn  according 
to  II g.. 

Consider  now  the  subproblem  of  using  the  sensor 
readings  from  a  sensor  at  location  If,  to  determine  if 
the  CBRN  release  corresponds  to  scenario  Li  or  Lj.  The 
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Generalized  Likelihood  Ratio  Test  (GLRT)  compares  the 
normalized  generalized  log-likelihood  ratio 


Xijk{yk'n) 


llo  suPo<gn<i%(y*’n) 

n  °Ssup0.en.  p]?.(yk’n) 


to  a  threshold  A  and  declares  Li  to  be  the  release  scenario 
whenever 


Furthermore,  assuming  (8)  is  satisfied 

limsup-log aijk,n(0j)GLRT  <  -A, 

n — >oo  R 

for  ail  Oj  G  flj,  and 

limsup  -  log <  -  inf  72(Q|n£ .), 

ra->oo  Tl  J  QG^yfc 


yk,n  €  y^T  =  {yn\Xijk(yn)  >  A},  (6) 

and  Lj  otherwise.  Selection  of  the  value  for  A  in 
(6)  is  performed  through  evaluation  of  decision  error 
probabilities  and  is  discussed  later  in  this  section. 

The  GLRT  described  in  (6)  can  make  two  types  of 
error.  Namely,  the  decision  test  can  declare  Lj  to  be  the 
release  scenario  when  in  actuality  the  release  scenario 
was  Li  and  vice  versa.  The  probabilities  of  these  errors 
are  represented  as, 

oGLRT  (a  \  _  td k  r  fc,n  a  cyxGLRT] 

Pijk,n  \yi)  ~  ^Qi  [y  T  ^ ijk,n  J • 

We  will  refer  to  the  quantities 

limsup-  log 
n — ►oo  R 

limsup-  log  Pfj%*T{9i), 
n — >oo  R 

as  the  exponents  of  the  Type-I  and  Type-II  error  proba¬ 
bilities,  respectively. 

Certainly,  some  values  of  A  in  (6)  are  better  than 
others.  To  evaluate  the  GLRT’s  performance,  we  will 
use  a  generalized  version  of  Definition  1,  which  we  will 
call  the  generalized  Neyman-Pearson  criterion. 

Definition  2 

The  decision  rule  {S^ \jk,n }  is  optimal  if  it  satisfies 

limsup  -  log otfjk  n{6j)  <  -A,  M6j  G  flj  (7) 

n — >oo  R 


for  all  9i  G  f 2j. 

Now,  in  the  event  that  the  GLRT  is  not  asymptotically 
optimal  (i.e.,  condition  (8)  does  not  hold)  information  on 
the  error  probabilities  can  still  be  gleamed.  Define  the  set 

%jk  =  {Q|  inf 0j(za:j  72 (Q |n^ )  —  inf ^ 72(Q|Il£.)  < 

A}.  Note  yk’n  G  ^^k^  if  and  only  if  £  S>ijk. 

From  Theorem  V.l  it  follows: 

limsup-log 

n— »oo  R 

<-  inf  72(Q|Il£  ), 

Q  t@ijk  1 

limsup-log  (9i) 

n— >oo  R 

<-  inf  /2(Q|n0feJ, 

Q  E@ijk 

In  addition,  we  have  the  following  lemma,  proven  in 
Appendix  B. 

Lemma  V.3  If  the  GLRT  is  not  optimal,  the  constraint 
limsup  -  log <  -A,  M9j  G  FL3 

n — »oo  R 

still  holds. 

We  can  therefore  determine  an  asymptotic  bound  on 
the  exponent  of  the  Type-II  error  probability  by  the 
nonlinear  program 

Z ijk  (A  ,9i)=  minQ  72(Q|II$.) 

s.t.  <  A,  (11) 

Aq  <  0, 


(9) 

V6*j  G  f lj 

(10) 

Wi  G  f 


and  maximizes  —  limsup„^00  ^  log (df^kn{9i)  uni¬ 
formly  over  all  9i  G  Q,. 

The  following  theorem,  whose  proof  is  in  Appendix 
A,  establishes  the  optimality  of  the  GLRT  according  to 
Definition  2. 

Theorem  V.2  The  GLRT  with  a  threshold  A  is  optimal 
under  Definition  2  if  and  only  if 

inf  J2(Q|Il£  )  >  inf  J2(Q|Ilg  ),  V^eOi, 

(8) 

where  ,Aljk  =  {Q|  inf^.gn^  72(Q|II^)  <  A}  and 
& 'ijk  =  {Q|inf^.eOj  /2(Q|n^)-infejefli  72(Q|Il£.)  < 
A  <  inf^eo,  /2(Q|n^.)}. 


where  #(Q;  fi.,-)  =  minify  72(Q|n^.)  - 

min^gQj  72(Q|Il0.).  The  constraints  enforcing  Q 
to  be  shift  invariant  and  correspond  to  valid  transition 
probabilities  for  a  Markov  chain  are  linear  and  are 
represented  by  the  system  of  linear  inequalities  Aq  <  0, 
where  q  is  a  “vectorization”  of  Q. 

Using  (11),  Zijk( A)  =  min^gn,  Zijk( A,  9f)  is  the  ex¬ 
ponent  of  the  worst  case  Type-II  error  probability  of  the 
GLRT  (i.e.,  declaring  the  release  scenario  is  Lj  when  the 
actual  release  scenario  is  Lf)  when  we  use  information 
from  sensor  k.  Zijkf )  is  a  non-increasing,  non-negative 
function,  with  lim^oc  Zijk(\)  =  0.  Therefore,  there 
must  exist  X*jk  >  0  such  that  Zijk( X*jk)  =  X*jk,  which 
represents  the  exponent  of  equal  Type  I  and  Type  II  error 
probabilities  in  the  GLRT.  Since  there  is  no  guarantee 
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that  X*jk  =  \*ik  (i.e.,  information  from  a  sensor  at 
location  Bk  does  not  provide  symmetric  decision  making 
via  the  GLRT),  we  define 

dijk  =  max{A*ifc,  A*ife},  (12) 

and  this  represents  a  symmetric  error  probability  expo¬ 
nent.  In  particular,  if  X*]k  is  the  maximizer  in  (12)  we  let 
(i,j)  =  otherwise  we  let  (i,j)  =  ( j,i ).  In  either 

case,  we  use  the  GLRT  test  that  compares  2Qjfc(yfe>™) 
to  dijk ■  The  above  leads  to  the  following  proposition. 

Proposition  V.4  Suppose  the  sensor  at  Bk  uses  the 
GLRT  and  compares  Xjjk{yk'n)  to  dijk ■  Then,  the 
maximum  probability  of  error  satisfies 

hmsup  i  log  P^l  <  - d^k , 

n—>oo  iL 
(e) 

where  PLk  n  denotes  the  maximum  probability  of  error 
when  deciding,  via  the  GLRT,  between  release  scenarios 
Li  and  Lj  using  information  from  a  sensor  placed  at 
location  Bk- 

Thus,  for  every  sensor  location  in  3d  we  have  an  upper 
bound  on  the  GLRT’s  probability  of  error  when  deciding 
between  any  two  scenarios  in  2z? .  This  result  will  prove 
useful  when  it  comes  to  selecting  locations  to  place  a 
limited  number  of  sensors. 

A.  Fusing  information  from  multiple  sensors 

So  far  all  of  the  results  presented  for  the  source 
localization  problem  considered  the  binary  GLRT  in  (6) 
when  the  real  problem  is  to  select  a  location  that  appears 
within  the  scenarios  in  3£ .  One  approach  is  to  make  a 
series  of  decisions  of  the  form  (6)  until  a  single  release 
scenario  has  been  chosen. 

Assume,  without  loss  of  generality,  that  the  sensors  are 
placed  at  locations  B i, . . . ,  Bk-  Using  the  observations 
from  sensor  Bk *,  where  k*  =  argmaxfcc(i2fc,  a  decision 
is  made  between  scenarios  L  t  and  L2  via  the  GLRT 
with  threshold  d12fc*,  as  established  in  (6)  and  (12).  This 
scenario  is  then  compared  with  L:s  by  the  GLRT,  and  so 
on,  until  N  —  1  decisions  have  been  made  and  only  one 
scenario  from  Jz?  is  accepted.  The  location  associated 
with  this  selected  scenario  is  declared  the  CBRN  source 
location.  We  will  refer  to  this  localization  paradigm  as 
the  sequential  GLRT  approach. 

Alternatively,  localization  can  be  performed  using  all 
sensor  information  in  unison  by  selecting  the  location 
with  greatest  likelihood  in  a  single  decision.  We  will 
refer  to  this  decision  paradigm  as  the  centralized  max¬ 
imum  likelihood  approach.  Unfortunately,  the  amount 
of  information  needed  to  fully  fit  a  complete  joint 
probability  mass  function  for  all  sensor  observation 
combinations  makes  it  practically  impossible.  Instead, 


joint  log-likelihood  values  under  scenario  Lj  can  be 
approximated  via  (4)  and 

K 

P(LJ)  =  ^log  sup  P£(yfc'").  (13) 

k=i  G  --L 

Localization  is  then  performed  by  selecting  the  location 
associated  with  L*  £  2zf  which  maximizes  the  value  in 
(13). 

The  sequential  GLRT  and  centralized  maximum  like¬ 
lihood  localization  techniques  can  be  combined.  In  what 
we  will  call  the  hybrid  approach  to  localization,  a 
sequence  of  decisions  is  made,  as  in  the  sequential 
GLRT  approach,  but  each  decision  employs  information 
from  more  than  one  sensors.  Assume,  again,  without 
loss  of  generality,  that  sensors  are  placed  at  locations 
B\,...,Bk-  Using  observations  from  m  sensors  at 
locations  Bk * , . . . ,  /j /,  ^  corresponding  to  the  m  largest 
values  of  d\2k,  k  =  1, . . . ,  K,  scenario  L\  is  declared  if 

771 

P(P1)=Vlog  sup  Pg*{yk*’n) 

£ _ ^ 

771 

>P(L2)=J2  log  sup  Pefe*V*>n), 

t—1  $£^2 

and  L2  is  declared  otherwise.  The  chosen  scenario  from 
this  decision  is  then  compared  to  scenario  L:i  in  the  same 
manner.  This  process  continues  until  only  one  release 
scenario  remains.  In  the  extreme  case  in  which  only 
one  sensor  is  used  in  each  test  within  the  sequence,  this 
approach  is  equivalent  to  the  sequential  GLRT  with  every 
threshold  set  to  0.  In  the  other  extreme  case,  in  which 
every  deployed  sensor  is  used  in  each  test  within  the 
sequence,  this  approach  is  equivalent  to  the  centralized 
maximum  likelihood  approach. 

VI.  Sensor  placement 

Now  that  performance  guarantees  have  been  estab¬ 
lished  for  each  sensor  location  in  33,  the  question  of 
where  to  place  TC  <  M  sensors  can  be  addressed. 
Using  the  bound  in  Proposition  V.4,  we  can  select  the 
K  elements  of  33  such  that  the  worst  case  error  when 
differentiating  between  any  two  release  scenarios  in  2z?  is 
minimized.  This  is  done,  as  in  [4],  via  the  mixed  integer 
linear  program  (MILP) 

max  e 

,  I/ 

/  vU— i  xk  —  K) 

Efc=l  Vijk  1?  f  j  1?  •  •  *  ?  t  X  j, 

Dijk  <  xk,  Vi,j,i  <j,k=  1, . . . ,  M,  (14) 

e  1  dijkUijki  i  <  j , 

Vijk  >  o,  Vi,j,k,i  <  j, 
xk  €  {0, 1},  Vfc. 

In  (14),  the  Xk  are  the  indicator  functions  for  place¬ 
ment  of  a  sensor  at  location  Bk .  Since  only  K  sensors 
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are  available,  the  number  of  positive  indicator  variables 
is  limited  in  the  constraint  ^2k=1  xk  —  K.  The  constraint 
e  <  Yjk=xdijk1hjk,  where  dijk  is  defined  in  (12), 
ensures  that  for  every  release  scenario  pair  in  Jif,  there 
exists  a  d^k  at  least  as  large  as  e.  Although  shown  to 
be  NP-hard,  (14)  can  be  solved  efficiently  for  large  sets 
SB  and  2zf  using  a  special  purpose  algorithm  [17]. 

Let  <3B  be  any  subset  of  SB  with  cardinality  I\ .  Further 
let  x('?^)  denote  the  vector  of  decision  variables  where 
Xki^S)  denotes  the  presence  of  Bk  in  W .  Define 

GfV)  =  min  max  dak- 
i<j  {fc|xfe(30=l} 

This  term  simply  represents  the  worst  case  pair-wise 
error  probability  exponent  present  in  the  placement  solu¬ 
tion  .  Allowing  x*,  y*,  and  e*  to  represent  the  optimal 
decision  variables  in  (14),  we  have  the  following  result 
from  [17],  but  now  with  dijk  s  corresponding  to  the 
Markovian  rather  than  the  iid  case. 

Proposition  VI.  1  For  any  sensor  placement  SJ  we  have 
e *  >  e(&). 

Moreover, 

e*  =  min  max  duk  (15) 

i<j  {k\xl~l} 

and  the  optimal  solution  satisfies 

*  f  L  if  k  =  argmaxWx,=1}dijk, 

dl-,k  |  0,  otherwise, 

for  all  i,  j,  i  <  j,  and  all  k,  where  at  most  one  y*-k  is 
set  to  1  for  a  given  (i,j)  pair. 

For  the  sequential  GLRT  decision  making  scheme,  a 
performance  guarantee  for  sensor  placement  determined 
via  solution  to  (14)  is  found  in  the  following  proposition. 
The  proof  is  provided  in  Appendix  C. 

Proposition  VI.2  Let  x*,  y* ,  and  e*  be  an  optimal 
solution  to  the  MILP  in  (14).  Place  sensors  at  loca¬ 
tions  Bk  such  that  xk  =  1  and  for  every  (Li,Lj) 
scenario  pair,  select  one  sensor  with  index  k*j  such  that 
yijk »  =  1.  Then,  the  worst  case  probability  of  error  for 
the  sequential  GLRT  decision  rule,  pjf^,opt,  satisfies 

limsup <  -e*. 

n — »oo  17 


VII.  Numerical  experimentation 

To  demonstrate  the  performance  of  the  sensor  place¬ 
ment  approach  of  Section  VI  and  the  two-stage  CBRN 
release  detection  and  localization  of  Section  III,  we 
simulated  several  release  scenarios  in  two  illustrative 


environments.  In  the  first  environment  -  the  toy  environ¬ 
ment  -  CBRN  releases  were  simulated  in  the  Quick  Ur¬ 
ban  &  Industrial  Complex  (QUIC)  Dispersion  Modeling 
System  [5]  developed  at  the  Los  Alamos  National  Lab¬ 
oratory.  QUIC  first  solves  the  fluid  dynamics  problem 
of  determining  local  wind  eddies  throughout  a  modeled 
three-dimensional,  outdoor  setting  using  the  methods  of 
Rockle  [18].  Using  the  fluid  flow  solution,  QUIC  simu¬ 
lates  the  travel  of  CBRN  particulates  via  a  Lagrangian 
random  walk.  Previously,  the  QUIC  codes  have  been 
tested  and  validated  for  real-world  situations  [5],  A  three- 
dimensional  illustration  of  the  toy  environment  some 
time  after  a  QUIC  simulated  point  release  is  found  in 
Figure  1. 

Additionally,  CBRN  releases  were  simulated  using  a 
simulator  we  developed  based  on  the  Lattice  Boltzmann 
Method  ( LBM ).  LBM  evolved  from  the  numerical  fluid 
modeling  technique  Lattice  Gas  Automata  (LGA),  in 
which  parcels  of  air  adhere  to  microscopic  laws  which 
dictate  their  movement.  Macroscopic  values  of  flow  ve¬ 
locities  and  densities  are  then  derived  by  the  underlying 
microscopic  properties  propagated  by  the  algorithm  [19]. 
Unfortunately,  LGA  often  falls  victim  to  instability  in 
the  face  of  statistical  noise  [20],  LBM  extends  LGA 
by  considering  air  parcel  movement  more  notionally  by 
modeling  microscopic  air  parcel  velocities  as  distribu¬ 
tions  in  the  Lattice  Boltzmann  Equation  (LBE).  It  has 
been  shown  that  under  reasonable  starting  conditions 
LBM  provides  accurate  approximations  to  fluid  flows. 
Further,  the  macroscopic  Navier-Stokes  equations  can  be 
recovered  from  the  microscopic  LBE  [6]. 


Fig.  1.  Toy  environment  QUIC  representation.  The  dots  represent 
dispersed  particulates  some  time  after  a  point  release.  The  release 
location  is  between  the  cylinder  and  the  top-most  block  in  the  figure. 

Wind  is  the  main  driver  of  particulate  transportation. 
As  such,  fluid  flow  accuracy  should  be  of  prime  concern 
when  selecting  a  dispersion  simulator.  This  becomes 
even  more  important  if  the  environment  of  interest  con¬ 
tains,  as  in  the  scenarios  considered  here,  obstructions 
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or  irregular  geometries.  These  characteristics  produce 
irregularities  in  the  environment’s  wind  field.  Both  QUIC 
and  LBM  accept  nominal  wind  conditions  which  are 
translated  into  varying  wind  eddies  distributed  through¬ 
out  the  environment. 

It  should  be  noted  that  Markov  chains  have  been 
shown  to  be  good  approximations  of  Lagrangian  ran¬ 
dom  walk  based  dispersion  models  [21],  such  as  the 
dispersion  modeling  used  by  QUIC.  One  should  expect 
to  see  our  localization  approach  to  perform  well  on 
QUIC  generated  data.  Our  LBM  model,  however,  is  not 
susceptible  to  the  same  ailment. 

We  stress  that  these  two  independently  developed 
simulators  make  none  of  the  distributional  assumptions 
(Markovian  particulate  concentrations)  our  localization 
methods  postulate.  They  simply  consider  a  point  release 
and  simulate  the  physics  of  dispersion  in  the  given 
environment  and  under  a  given  weather  pattern.  Hence, 
the  results  we  produce  assess  not  only  the  performance 
of  the  methods  but  also  the  validity  of  our  assumptions. 

In  the  second  environment  -the  city  grid  environment- 
CBRN  releases  were  also  simulated  in  QUIC  and  LBM. 
This  environment  consists  of  geometries  more  typical 
of  dense  urban  areas.  It  was  within  this  environment 
that  we  compared  the  agent  concentration  profiles  of 
the  two  different  dispersion  simulators  at  a  sampler 
located  down  wind  of  the  point  release.  This  compar¬ 
ison  appears  in  Figure  2.  These  models  have  different 
discretizations  of  the  three-dimensional  model  and  hence 
produce  concentration  values  that  differ  in  scale.  This 
has  been  accommodated  in  Figure  2  by  reporting  the 
percentage  of  all  observed  concentrations  reported  by 
a  single  sensor  downwind  from  a  release.  The  LBM 
simulator  produces  a  much  smoother  agent  concentration 
evolution  than  the  QUIC  simulator.  Thus,  noise  within 
the  evaluation  of  the  proposed  methodologies  when  the 
LBM  data  are  used  is  primarily  due  to  our  sensor  false 
alarm  model  rather  than  the  dispersion  model.  Further 
dissimilarities  between  the  two  concentration  profiles 
also  stem  from  differences  in  the  way  surface  conditions 
are  treated  in  the  two  simulated  environments.  QUIC 
employs  boundary  conditions  to  reflect  interaction  with 
a  concrete  surface.  The  LBM  representation,  on  the  other 
hand,  uses  only  the  so-called  “bounce  back”  boundary 
conditions  [6], 

A.  Toy  Environment 

To  demonstrate  the  performance  of  the  proposed 
techniques,  first  a  simple  three-dimensional  environment 
with  several  large  obstructions  (i.e.,  buildings  in  an 
urban  setting)  was  modeled  in  QUIC  and  LBM.  CBRN 
point  releases  originating  from  three  release  locations 
were  simulated  with  wind  originating  from  NNW,  N,  or 
NNE.  225  potential  sensor  locations  are  considered  on  an 


Sensor  1 2  Concentration  Profile 


Time  Step 


Fig.  2.  Evolution  of  concentration  at  a  sensor  downwind  from  the 
source  of  the  CBRN  event. 


evenly  spaced  grid  near  ground  level.  Figure  3  illustrates 
a  cross  section  of  the  environment  with  the  locations  of 
the  obstructions  (larger  shapes)  and  releases  (denoted  by 
an  “x”). 


0' - 1 - 1 - 1 - 1 - 1 - 1 

0  5  10  15  20  25  30 

Meters 

Fig.  3.  Toy  environment  with  optimal  sensor  placement  solution. 

The  probability  laws  were  determined  empirically  via 
Monte  Carlo  simulation  (using  (1)  and  (5))  for  each 
combination  of  wind  direction,  release  location,  and 
sensor  location.  First,  sensor  concentration  observations 
are  encoded  into  binary  outputs  indicating  the  presence 
of  non-nil  concentration.  Considering  the  binary  nature 
of  the  sensor  readings,  sensor  measurement  noise  is 
introduced  via  Bernoulli  trials  which  result  in  flipping 
sensor  observations.  The  resulting  noisy  sensor  readings 
provide  a  wealth  of  observations  from  which  probability 
laws  can  be  established,  allowing  for  solution  of  (11) 
and  the  determination  of  the  dljk  in  (12).  For  evalu- 
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ation  of  the  release  detection  and  source  localization 
methodologies,  sensors  were  sampled  once  per  second 
for  50  seconds,  resulting  in  an  observation  window  of 
50  binary  indicators  of  the  presence  of  a  CBRN  agent. 
In  the  case  of  both  the  QUIC  and  LBM  data,  by  the 
time  six  sensors  are  placed  within  the  environment  there 
is  no  further  improvement  in  the  asymptotic  worst  case 
localization  error.  The  solution  depicted  for  the  case  in 
which  QUIC  data  were  used  is  shown  in  Figure  3  as 
the  small  circles.  The  small  squares  correspond  to  the 
optimal  sensor  placement  determined  using  LBM  data. 

A  test  set  is  constructed  for  each  release  location  by 
first  selecting  a  wind  direction  according  to  a  “wind 
rose”  which  describes  the  likelihood  of  each  wind  di¬ 
rection.  Then,  using  the  simulations  from  which  the 
probability  laws  were  derived,  test  sensor  observations 
are  generated. 

1 )  Release  detection  evaluation:  For  each  scenario 
in  the  test  set,  a  release  was  declared  whenever  at 
least  one  sensor’s  observations  resulted  in  declaring  a 
CBRN  attack  according  to  the  test  in  (2).  To  gain  an 
understanding  on  the  CBRN  detector’s  performance,  this 
process  was  conducted  for  many  different  values  of  :/ 
to  facilitate  the  construction  of  the  Receiver  Operating 
Characteristic  (ROC)  curve.  When  a  reasonable  sensor 
error  rate  of  1/8  was  used  in  both  the  case  in  which 
QUIC  data  were  used  in  the  evaluation  and  when  LBM 
data  were  used,  the  detector  attained  perfect  detection 
performance,  with  100%  probability  of  detection  and 
0%  probability  of  false  alarm.  To  better  understand  how 
robust  our  detection  methodology  is  to  sensor  noise,  we 
repeated  the  experiment  with  different  sensor  error  rates. 
Figure  4  shows  the  probability  of  detection  achieved 
when  0%  probability  of  false  alarm  is  demonstrated 
as  a  function  of  sensor  error  rate.  As  shown,  there  is 
a  precipitous  decline  in  detection  performance  by  the 
time  roughly  one  third  of  the  sensor  measurements  are 
erroneous. 

2)  Source  localization  evaluation:  The  ability  of  all 
three  localization  methodologies  we  outlined  in  Sec¬ 
tion  V  was  also  evaluated  on  each  test  scenario.  Table  I 
displays  the  accuracy  of  each  approach. 

TABLE  I 

Accuracy  of  localization  methodologies  in  the  toy 

ENVIRONMENT 


QUIC  Data 

LBM  Data 

Sequential 

99.7% 

98.3% 

Centralized 

100% 

99.7% 

Hybrid 

100% 

99.7% 

When  the  hybrid  localization  approach  is  used,  gains 
in  localization  accuracy  begin  to  disappear  after  a  certain 
number  of  sensors  are  used  in  each  decision  within  the 
sequence  of  maximum  likelihood  tests.  This  accuracy 


Probability  of  Detection  by  Two  State  Sensor  Error  Rate 


Fig.  4.  Probability  of  detection  as  the  error  rate  for  the  two  state 
CBRN  sensor  grows.  Reported  values  correspond  to  zero  observed 
false  alarms. 


is  equal  to  the  accuracy  attained  by  the  centralized 
maximum  likelihood  approach.  When  QUIC  data  were 
used  for  evaluation,  perfect  performance  was  observed 
when  only  one  sensor  was  used  in  each  decision  in  the 
sequence.  This  actually  corrected  the  one  mistake  found 
when  the  sequential  GLRT  localization  technique  was 
used.  When  LBM  data  were  used  for  evaluation,  it  was 
not  until  three  sensors  were  employed  for  each  decision 
within  the  sequence  of  maximum  likelihood  tests  that  the 
same  performance  as  the  maximum  likelihood  approach 
was  observed. 

B.  City  Grid  Environment 

To  further  evaluate  the  proposed  CBRN  detector  and 
locator,  a  three-dimensional  environment  inspired  by  a 
dense  urban  grid  pattern  was  modeled  in  QUIC  as  well  as 
using  LBM.  The  environment  consists  of  a  city  grid,  four 
blocks-by-four  blocks.  Each  block  is  100  meters-by-100 
meters  with  10  meter-wide  throughways.  Each  block’s 
height  is  drawn  randomly  from  the  uniform  distribution 
ranging  from  20  to  60  meters.  Sensors  are  allowed  to 
be  placed  at  any  intersection  and  five  intersections  are 
considered  as  potential  release  locations.  The  shape  of 
the  grid,  as  well  as  the  location  of  simulated  releases,  is 
shown  in  Figure  5. 

CBRN  event  scenarios  considered  are  the  40  unique 
combinations  of  five  different  release  locations  within  the 
grid  with  wind  blowing  at  1  m/s  or  5  m/s  and  originating 
from  the  four  cardinal  directions.  The  two  selected  wind 
speeds  represent  nominal  wind  conditions,  but  do  not 
translate  into  a  uniform  wind  field.  As  the  simulation 
proceeds,  local  wind  eddies  vary  according  to  conditions 
incited  by  the  presence  of  buildings  within  the  grid. 
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Sensor  Location:  QUIC  data 
Sensor  Location:  LBM  data 
Release  location 


A 


Fig.  5.  City  environment  with  optimal  sensor  placement  solution. 


As  in  Subsection  VII-A,  the  probability  laws  ac¬ 
cording  to  (1)  and  (5)  for  each  combination  of  CBRN 
scenario  and  sensor  location  are  determined  by  Monte 
Carlo  simulation  using  data  generated  from  both  QUIC 
and  LBM  simulations.  This  time,  however,  not  only  is 
the  binary  sensor  model  with  Bernoulli  false  alarms 
considered,  but  a  real  valued  sensor  model  with  additive 
white  noise  is  considered  as  well.  These  sensors  are 
modeled  as 

C  =  C  +  N(0,ae),  (16) 

where  C  is  the  actual  concentration  of  particulate  present 
at  the  sensor’s  location,  7V(0,<re)  denotes  a  normally 
distributed  random  variable  with  mean  0  and  standard 
deviation  ae,  and  C  denotes  the  sensor’s  reported  con¬ 
centration  observation.  Values  of  C  were  then  translated 
into  a  discretized  representation  of  five  concentration 
states.  First  a  value  of  <jt  which  produces  an  equivalent 
false  alarm  as  the  binary  sensor  model  was  used.  To 
illustrate  the  power  of  increasing  the  number  of  concen¬ 
tration  states  used  in  the  analysis,  the  normal  random 
variate  in  (16)  was  replaced  with  iV(0, 10cre). 

For  the  binary  sensor  case,  if  data  generated  via 
LBM  were  used  for  analysis,  little  improvement  in  worst 
case  error  probability  is  observed  after  six  sensors  are 
deployed.  However,  when  data  generated  via  QUIC  were 
used  for  analysis,  it  is  not  until  16  sensors  are  deployed 
that  no  further  improvement  is  made  in  the  bound  on  the 
worst  case  error  probability. 

For  the  real-valued  sensor  case  with  mild  noise,  if 
data  generated  via  LBM  were  used  for  analysis,  once 
five  sensors  are  placed  in  the  environment  no  further 
improvement  in  worst  case  error  probability  is  observed. 
These  five  sensor  locations  correspond  to  the  five  release 
locations  in  the  environment.  When  QUIC  data  were 


employed  for  analysis,  no  further  improvement  in  the 
worst  case  error  probability  is  observed  after  13  sensors 
have  been  deployed.  The  optimal  placements  of  real¬ 
valued  sensors  with  mild  noise  are  presented  in  Figure  5. 

In  the  case  of  real  valued  sensor  models  with  large 
noise,  analysis  utilizing  LBM  data  shows  no  further 
improvement  after  eight  sensors  are  deployed.  When 
QUIC  data  is  used,  no  further  performance  is  guaranteed 
after  the  same  13  sensors  are  deployed.  This  solution  is 
the  same  as  the  one  found  in  the  mild  noise  case. 

A  test  set  for  the  city  grid  was  constructed  for  each 
of  the  release  locations  by  first  selecting  a  wind  direc¬ 
tion  (all  four  cardinal  directions  are  chosen  with  equal 
likelihood)  and  then  selecting  a  wind  speed  (both  wind 
speeds  are  selected  with  equal  likelihood).  The  test  set 
sensor  readings  were  generated  in  the  same  way  as  the 
sensor  readings  used  to  determine  the  probability  laws. 

1 )  Release  detection  evaluation:  Following  the  proce¬ 
dure  in  Subsection  VII-A1,  for  each  scenario  in  the  test 
set,  a  release  was  considered  detected  whenever  at  least 
one  sensor’s  observations  lead  to  the  decision  of  “attack” 
via  the  test  in  (2).  The  resulting  ROC  curves  are  depicted 
in  Figure  6.  Whether  QUIC  or  LBM  simulations  were 
used  to  generate  scenario  data,  the  detector  attains  nearly 
perfect  performance.  This  holds  true  regardless  of  which 
sensor  model  is  used. 
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Fig.  6.  ROC  curve  of  CBRN  detector  in  the  city  grid  environment. 


2)  Source  localization  evaluation:  The  proposed 
CBRN  localization  techniques  performed  well  on  both 
the  QUIC  and  LBM  generated  data  sets,  even  when  real 
valued  sensors  with  large  additive  noise  are  used.  The 
performances  observed  are  listed  in  Table  II. 

While  all  three  methodologies  perform  phenomenally 
in  both  simulated  environments,  the  hybrid  and  central¬ 
ized  maximum  likelihood  techniques  are  the  apparent 
champions,  both  consistently  outperforming  the  sequen¬ 
tial  localization  technique.  Yet,  the  sequential  technique 
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has  the  advantage  of  being  fully  distributed  which  is 
desirable  in  implementations  involving  a  large  number 
of  sensors. 

In  some  instances,  the  hybrid  approach  outperformed 
the  centralized  approach,  as  noted  in  the  binary  sensor 
and  large  noise  real  valued  sensor  scenarios  using  LBM 
data.  In  these  cases,  reducing  the  number  of  sensors  used 
in  a  maximum  likelihood  decision  prevented  a  single 
sensor’s  observations  skewing  a  computed  probability 
into  an  incorrect  localization. 


C.  City  Grid  Environment  with  Sensor  Constraints 

In  Subsection  VII-B,  the  optimal  solution  to  (14)  was 
evaluated  when  K  was  determined  by  performance  gain 
analysis.  However,  it  may  be  the  case  that  a  limited 
number  of  sensors  are  available  for  deployment.  To 
assess  the  performance  of  the  proposed  CBRN  detector 
and  source  locators  in  a  sensor  constrained  scenario,  we 
repeated  the  analysis  in  Subsection  VII-B,  this  time  with 
K  =  5,  indicating  that  five  sensors  are  available  for 
deployment.  The  optimal  placements  of  binary  sensors 
under  this  constraint  can  be  found  in  Figure  7  using 
the  QUIC  and  LBM  generated  data  and  binary  sensor 
models. 


Fig.  7.  Optimal  placement  of  5  sensors  in  the  city  environment. 

When  the  number  of  sensors  deployed  decreases,  the 
number  of  potentially  invisible  releases  (i.e.,  a  release 
that  is  not  detected  by  a  single  sensor)  increases.  An 
intuitive  approach  to  address  this  concern  is  to  span  the 
deployed  sensors  out  as  much  as  possible.  It  appears,  in 
Figure  7,  that  the  optimal  placement  solutions  adhere  to 
this  philosophy.  Evaluations  using  the  real  valued  sensor 
models  produced  optimal  sensor  placements  correspond¬ 
ing  to  the  five  simulated  release  locations. 


1 )  Release  detection  evaluation:  There  is  a  significant 
decrease  in  performance  due  to  the  increased  presence 
of  invisible  releases  in  the  test  sets,  as  illustrated  in 
Figure  8.  However,  this  degradation  is  confined  solely  to 
analysis  which  uses  a  binary  sensor  model  and  CBRN 
dispersion  modeled  using  QUIC. 

ROC:  CBRN  Release  Detection 


Fig.  8.  ROC  curve  of  CBRN  detector  in  the  city  grid  environment 
for  the  case  of  using  only  five  sensors. 

2)  Source  localization  evaluation:  Localization  ac¬ 
curacy  for  all  three  approaches  is  listed  in  Table  III. 
We  actually  see  an  increase  in  performance  when  the 
number  of  sensors  are  limited  in  the  LBM  scenarios 
with  real  valued  sensors  and  large  additive  noise.  This  is 
due  to  removing  sensor  information  that  previously  led 
to  erroneous  localizations,  a  phenomenon  similar  to  the 
reason  we  observed  better  performance  from  the  hybrid 
localization  technique  over  the  centralized  approach. 

An  apparent  follow  up  question  to  the  localization 
performance  analysis  is,  how  well  does  the  optimal 
solution  to  the  sensor  placement  MILP  compare  to 
other,  non-optimal  solutions?  To  answer  this  question, 
placements  in  which  the  sensor  locations  are  chosen 
uniformly  randomly  were  evaluated  on  the  same  test  set 
as  the  optimal  placement  solution.  The  percentage  of  test 
cases  that  were  correctly  localized  appear  in  Figure  9. 
Clearly,  the  optimal  placement  solution  outperforms  any 
of  the  randomly  chosen  sensor  placements. 

VIII.  Conclusions 

We  have  presented  a  new  two-stage  methodology  to 
the  related  problems  of  detecting  CBRN  events  and  then 
determining  the  location  of  their  origin.  Our  approach 
bypasses  many  of  the  problems  “inverse  dispersion” 
approaches  face  in  an  inhomogeneous  and  uncertain 
urban  setting. 
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TABLE  II 

Accuracy  of  localization  methodologies  in  the  city  grid  environment 


QUIC  Data 
Binary  sensors 

QUIC  Data 

Real  valued  sensors 
Mild  noise 

QUIC  Data 

Real  valued  sensors 
Large  noise 

LBM  Data 
Binary  Sensors 

LBM  Data 

Real  valued  sensors 
Mild  noise 

LBM  Data 

Real  valued  sensors 
Large  noise 

Sequential 

98.42% 

93.50% 

94.66% 

98.98% 

94.18% 

83.90% 

Centralized 

99.98% 

100% 

99.98% 

99.82% 

100% 

99.48% 

Hybrid 

99.98% 

100% 

99.32% 

99.84% 

100% 

99.60% 

TABLE  III 

Accuracy  of  localization  methodologies  in  the  city  grid  environment  when  only  five  sensors  are  available  for 

deployment 


QUIC  Data 
Binary  sensors 

QUIC  Data 

Real  valued  sensors 
Mild  noise 

QUIC  Data 

Real  valued  sensors 
Large  noise 

LBM  Data 
Binary  Sensors 

LBM  Data 

Real  valued  sensors 
Mild  noise 

LBM  Data 

Real  valued  sensors 
Large  noise 

Sequential 

90.76% 

93.40% 

93.90% 

98.72% 

94.18% 

86.48% 

Centralized 

90.02% 

96.74% 

96.62% 

98.66% 

100% 

99.52% 

Hybrid 

90.02% 

97.88% 

99.32% 

97.40% 

100% 

99.60% 

Random  Sensor  Placement  Performance 


0.9 


O.f 


0 


0.7 


- QUIC  data,  optimal  placement,  sequential  GLRT 

0.3  - - QUIC  data,  optimal  placement,  maximum  likelihood 

LBM  data,  optimal  placement,  sequential  GLRT 

- LBM  data,  optimal  placement,  maximum  likelihood 

g  2  O  QUIC  data,  sequential  GLRT 

*  QUIC  data,  maximum  likelihood 
□  LBM  data,  sequential  GLRT 
LBM  data,  maximum  likelihood 

0.1 1  '  '  '  '  '  — i - 1 - c_ 

0  2  4  6  8  10  12  14  16 

Random  Solution  Number 


8  8 

o 


decision  making  process,  yet  it  offers  obvious  implemen¬ 
tation  advantages  in  large  deployments.  The  centralized 
maximum  likelihood  localization  approach  remedies  this 
by  employing  all  deployed  sensor  information  in  a  single 
decision,  however,  the  hybrid  approach  committed,  in  to¬ 
tal,  slightly  fewer  localization  errors  in  our  experiments. 

Appendix  A 
Proof  of  Theorem  V.  2 

Our  approach  to  show  optimality  of  the  GLRT  will  be 
to  show  that  it  is  a  subset  of  decision  rules  which  satisfy 
the  criterion  of  Definition  2.  To  that  end  we  first  present 
a  lemma  similar  to  a  result  first  produced  by  Hoeffding 
[16]. 


Fig.  9.  Evaluation  of  20  randomly  chosen  placements  of  5  sensors. 


To  detect  attacks,  relative  entropy  is  used  to  compare 
current  sensor  observations  to  known  nominal  probabil¬ 
ity  laws  derived  from  simulations.  Once  a  CBRN  event 
is  detected,  a  localization  procedure  is  used  to  deter¬ 
mine  where  the  CBRN  event  started.  Our  localization 
techniques  include  a  sequence  of  GLRT  decisions  using 
information  from  a  single  sensor,  a  maximum  likelihood 
selector  in  which  information  from  all  deployed  sensors 
are  used,  and  a  hybrid  of  these  two  approaches  consisting 
of  a  sequence  of  maximum  likelihood  decisions  using 
information  from  more  than  one  sensor. 

The  proposed  CBRN  detection  approach  demonstrated 
promise  within  a  numerical  evaluation  with  data  gen¬ 
erated  in  two  different  dispersion  simulators.  Our  ap¬ 
proaches  for  source  localization  performed  well  also. 
The  sequential  GLRT  localization  methodology  is  vul¬ 
nerable  to  committing  a  localization  error  early  on  in  the 


Lemma  A.l  The  sequence  of  tests 

Tjk,n  =  {yfe,nl  |ng )  >  a}  (I?) 

tfj  tiij 

is  optimal  under  Definition  2. 

Proof:  Consider  the  sequence  of  tests 

=  (yfc,"iwf  "|n£.)  >  A}. 

By  Lemma  1  of  [22],  for  all  9j  £  i  l:i ,  the  constraint 

limsup  -  log &ijk,v.(@j)  <  -A 

n — >oo  fl 

is  satisfied  and  the  value 

-  limsup  -  log 0ijk,n(Oi) 

n — >oo  ^ 

is  maximized  uniformly  over  all  6i  £  f where 
ctijk,n{0j)  and  (3ijk,n{9i)  indicate  the  Type-I  and  Type- 
II  errors  of  test  JT’ijk  n(& j) ■  Thus,  the  constraint  (7)  is 
satisfied  for  test  n,  since  yfc’n  £  y*-k  n  if  and  only 
if 

U 
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Further,  by  the  same  reasoning  and  using  the 
union  bound,  it  must  be  the  case  that  the  quantity 

-limsup^^  L  log /3,*j7Cj „(&,;)  is  maximized  uniformly 
over  all  9i  £  11,.  ■ 

Now,  allow  yk’n  £  y ’ijknT ■  Then  we  have, 

A  <  -  log  sup  Pg  (yk,n)  -  -  log  sup  Pg  (yfc,n) 
n  e^cii  n  e^Qj  J 

=  -  -  mfnh(£fn\Ilkei) 

--log  sup  Pg  (yk,n) 
n  0,-efL 

<  -  H(Syfn)  -  -  log  sup  Pg  (yk,n) 
n 

=  inf  I2(Sfn  \Ukej), 

tfj  fci  lj 

where  H(-)  is  the  empirical  entropy,  defined  as 

|e|  |£| 

(crv,au)  log  S1^  {au\av). 

U— 1 V=1 

The  first  and  second  equalities  stem  from  the  definitions 
of  H  and  I2 .  The  second  inequality  is  the  result  of 
the  non-negativity  of  I2 .  Given  (17),  it  follows  that 
^ijknT  —  & Tjk  n  and  Lemma  A.l  establishes  that  the 
constraint  (7)  is  satisfied  for  the  GLRT. 

For  any  9i  £  O,  note 

=PeAyk'n  i  Sffi* T } 

=PeAyk’n  i  ^k,n) 

+  [yk’n  £  ^ijk,ro  yk’n  £  ^ijk,nT]- 
If  yk,n  t  y&L»T, 

A  >  -  log  sup  Pg  (yk’n)  -  -  log  sup  Pgk( yk’n) 
n  ejgfjj  n  gjGnj  J 

=  -H{stn)-e  infy2(<ff"|n^) 

+  H(£y2k’n)+  inf  /2(<ff"  |n9j) 

tfj  Giij 

=  inf  I2{gf'n\Ue)  -  inf  I2(sfU  \H6i). 

Since  yk’n  £  y ’*jkn  implies  A  < 

mfg^h^K.),  when  yk’n  £  y^n  and 
yfc,n  ^  yGLRT  we  jjave  gy  ’  g  %jk.  By  Theorem 
V.l,  for  all  el  £  Slu 

-  limsup  -  logP£  [yfc’n  €  y*jk}n,  yk,n  i  y’Ajk,nT] 

>  inf  J2(  Q|H*). 

CJG  t &ijfc 

Therefore, 

-  limsup  -  log/3 AjktnT(Vi) 

n—>  00  Tt 


=  min  I  -  lim sup  -  log  (3*jkn (0»), 

n — kxj  ri 

-  limsup  -  logP^.  [yk,n  £  yt*jk,n>  yk'n  i  y'ijk,nT}\ 

n — >oo  Tl  J 

>min{  inf  /2(Q|n£)  inf  J2(Q|nfc)} 

=  inf  /2(Q|Hg  ). 

Q  £y*ijk 

It  follows  that  the  Type-II  error  probability  decreases 
exponentially  to  zero  faster  than  the  Type-II  error  prob¬ 
ability  of  test  y*jk  n.  Lemma  A.l  has  already  estab¬ 
lished  that  y*jk  n  is  optimal  under  Definition  2.  We 
can  therefore  conclude  that  the  GLRT  is  optimal  under 
Definition  2. 

Appendix  B 
Proof  of  Lemma  V.3 

If  the  GLRT  is  not  optimal,  for  all  Q  <f  $>ijk, 
inf^eo,  /2(Q|II^)  -  inf^eQ.  J2(Q|II^)  >  A.  Equiv¬ 
alently, 

inf  _  J2(Q|Ilg  )  >  A  +  inf  /2(Q|n*  )  >  A, 

tf  j  t  u  u  j  tf  ^  l  % 

where  the  last  inequality  follows  from  the  non-negativity 
of  J2.  Therefore,  using  Lemma  A.l,  it  follows  that 

limsup  —  log otfjkRT(9j)  <  —A,  \/9j  £  %. 

n — »oo  Tl 

Appendix  C 

Proof  of  Proposition  VI. 2 

The  sensor  with  index  k*,;  will  use  the  GLRT  and 
achieve  a  maximum  probability  of  error  with  exponent 
no  smaller  than  dljk*  .  For  every  i  and  j  7^  *  define 
En(i,j)  as  the  event  that  the  GLRT  employed  by  the 
sensor  at  Bk *  will  decide  L,  when  L,  is  true.  For  all 

ij  J 

5n  >  0  and  large  enough  n  we  have 
P 8i  [error]  <Pfti[U j7HEn(i,j)\ 

<Y^e~n{diik»+Sn)  <  (N  -  l)e-n{e*+Sn). 

The  second  inequality  above  is  due  to  Proposition  V.4 
and  the  last  inequality  above  is  due  to  (15).  Since  the 
bound  above  holds  for  all  i  we  obtain  the  desired  result. 

References 

[1]  G.  Schweitzer  and  E.  Sharber,  A.C.,  Countering  urban 
terrorism  in  Russia  and  the  United  States:  proceedings  of  a 
workshop.  National  Academies  Press,  2006.  [Online].  Available: 
http://books. google. com/books?id=qlW4eDiSvpgC 

[2]  T.  Cover  and  J.  Thomas,  Elements  of  Information  Theory.  Wiley, 
New  York,  1991. 

[3]  A.  Dembo  and  O.  Zeitouni,  Large  Deviations  Techniques  and 
Applications,  2nd  ed.  NY:  Springer- Verlag,  1998. 

[4]  I.  C.  Paschalidis  and  D.  Guo,  “Robust  and  distributed  stochastic 
localization  in  sensor  networks:  Theory  and  experimental  results,” 
ACM  Trans.  Sensor  Networks,  vol.  5,  no.  4,  pp.  34:1-34:22,  2009. 


LOCKE  AND  PASCHALIDIS:  MODEL-FREE  STOCHASTIC  LOCALIZATION  OF  CBRN  RELEASES 


15 


[5]  E.  R.  Paryak  and  M.  J.  Brown,  “QUIC  URB  v.1.1:  theory  and 
users  guide,”  Los  Alamos  National  Laboratory,  Tech.  Rep.,  2007. 

[6]  S.  Chen  and  G.  D.  Doolen,  “Lattice  Boltzmann  method  for  fluid 
flows,”  Annual  Review  of  Fluid  Mechanics,  vol.  30,  no.  1,  pp. 
329-364,  1998. 

[7]  A.  Atalla  and  A.  Jeremic,  “Estimation  of  boundary  properties  us¬ 
ing  stochastic  differential  equations,”  in  Electrical  and  Computer 
Engineering,  2009.  CCECE  ’ 09 .  Canadian  Conference  on,  2009, 
pp.  385-389. 

[8]  M.  Ortner  and  A.  Nehorai,  “A  sequential  detector  for  biochem¬ 
ical  release  in  realistic  environments,”  Signal  Processing,  IEEE 
Transactions  on,  vol.  55,  no.  8,  pp.  4173^1182,  2007. 

[9]  M.  Ortner,  A.  Nehorai,  and  A.  Jeremic,  “Biochemical  transport 
modeling  and  Bayesian  source  estimation  in  realistic  environ¬ 
ments,”  Signal  Processing,  IEEE  Transactions  on,  vol.  55,  no.  6, 
pp.  2520-2532,  2007. 

[10]  A.  Atalla  and  A.  Jeremic,  “Localization  of  chemical  sources  using 
stochastic  differential  equations,”  in  Acoustics,  Speech  and  Signal 
Processing  (ICASSP),  2008.  IEEE  International  Conference  on, 
2008,  pp.  2573-2576. 

[11]  M.  Ortner  and  A.  Nehorai,  “Biochemical  transport  modeling, 
estimation  and  detection  in  realistic  environments,”  in  Acoustics, 
Speech  and  Signal  Processing  (ICASSP),  2008.  IEEE  Interna¬ 
tional  Conference  on,  2008,  pp.  5169-5172. 

[12]  L.  Delle  Monache,  J.  K.  Lundquist,  B.  Kosovic,  G.  Johannesson, 
K.  M.  Dyer,  R.  D.  Aines,  F.  K.  Chow,  R.  D.  Belles,  W.  G.  Hanley, 
S.  C.  Larsen,  G.  A.  Loosmore,  J.  J.  Nitao,  G.  A.  Sugiyama,  and 
P.  J.  Vogt,  “Bayesian  inference  and  Markov  chain  Monte  Carlo 
sampling  to  reconstruct  a  contaminant  source  on  a  continental 
scale,”  Journal  of  Applied  Meteorology  and  Climatology,  vol.  47, 
no.  10,  pp.  2600-2613,  2008. 

[13]  A.  Keats,  E.  Yee,  and  F.  Lien,  “Bayesian  inference  for  source 
determination  with  applications  to  complex  urban  environment,” 
Atmospheric  Environment,  vol.  41,  pp.  465^479,  2007. 

[14]  M.  S.  Daskin,  Network  and  discrete  location:  models,  algorithms, 
and  applications,  ser.  Wiley-Interscience  series  in  discrete  math¬ 
ematics  and  optimization.  Wiley,  1995. 

[15]  G.  Hoblos,  M.  Staroswiecki,  and  A.  Aitouche,  “Optimal  design 
of  fault  tolerant  sensor  networks,”  in  International  Conference 
on  Control  Applications,  vol.  9,  no.  50,  2000,  pp.  467^172. 

[16]  W.  Hoeffding,  N.  Fisher,  and  P.  Sen,  The  collected  works  of 

Wassily  Hoeffding,  ser.  Springer  series  in  statistics:  Perspectives 
in  statistics.  Springer- Verlag,  1994.  [Online].  Available: 

http  ://books .  google .  com/books  7id=f7xLrFT  tKp  Y  C 

[17]  S.  Ray,  W.  Lai,  and  I.  C.  Paschalidis,  “Statistical  location 
detection  with  sensor  networks,”  Joint  special  issue  IEEE/ACM 
Trans.  Networking  and  IEEE  Trans.  Information  Theory,  vol.  52, 
no.  6,  pp.  2670-2683,  2006. 

[18]  Rockle,  “Bestimmung  der  stromungsverhlltnisse  im  bereich  kom- 
plexer  bebauungsstrukturen,”  Ph.D.  dissertation,  Technischen 
Hochschule,  Darmstadt,  Germany,  1990. 

[19]  U.  Frisch,  B.  Hasslacher,  and  Y.  Pomeau,  “Lattice-gas  automata 
for  the  Navier-Stokes  equation,”  Phys.  Rev.  Lett.,  vol.  56,  pp. 
1505-1508,  1986. 

[20]  P.  Lallemand  and  L.-S.  Luo,  “Theory  of  the  lattice  Boltzmann 
method:  Dispersion,  dissipation,  isotropy,  Galilean  invariance, 
and  stability,”  Phys.  Rev.  E,  vol.  61,  pp.  6546-6562,  2000. 

[21]  B.  Legg  and  M.  R.  Raupach,  “Markov-chain  simulation  of  par¬ 
ticle  dispersion  in  inhomogeneous  flows:  the  mean  drift  velocity 
induced  by  a  gradient  in  eulerian  velocity  variance,”  Boundary- 
Layer  Meteorology,  vol.  24,  pp.  3-13,  1982. 

[22]  O.  Zeitouni,  J.  Ziv,  and  N.  Merhav,  “When  is  the  generalized 
likelihood  ratio  test  optimal?”  Information  Theory,  IEEE  Trans¬ 
actions  on,  vol.  38,  no.  5,  pp.  1597-1602,  1992. 


PR.  Taylor  Locke  has  an  M.S.  in  Applied 
Mathematics  from  the  Rensselaer  Polytechnic 
Institute  (RPI)  and  a  Ph.D.  in  Systems  En¬ 
gineering  from  Boston  University  (BU).  He 
is  currently  at  the  MIT  Lincoln  Laboratory 
as  Technical  Staff  in  the  Informatics  and 
Decision  Support  Group. 

_ 


Ioannis  Ch.  Paschalidis  (M’96,  SM’06)  is  a 
Professor  and  Distinguished  Faculty  Fellow 
at  Boston  University  with  appointments  in 
the  Department  of  Electrical  and  Computer 
Engineering  and  the  Division  of  Systems 
Engineering.  He  is  a  Co-Director  of  the  Cen¬ 
ter  for  Information  and  Systems  Engineering 
(CISE).  He  completed  his  graduate  education 
at  the  Massachusetts  Institute  of  Technology 
(MIT)  receiving  an  M.S.  (1993)  and  a  Ph.D. 
(1996),  both  in  Electrical  Engineering  and 
Computer  Science.  In  September  1996  he  joined  Boston  University 
where  he  has  been  ever  since.  He  has  held  visiting  appointments 
with  MIT,  and  the  Columbia  University  Business  School.  His  current 
research  interests  lie  in  the  fields  of  systems  and  control,  networking, 
applied  probability,  optimization,  operations  research,  computational 
biology,  and  medical  informatics. 


